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The  self-consistent  electronic  structure  of  Si,  Ge,  and  zinc- 


blende  GaP,  GaAs ,  2nS  and  ZnSe 


determined  using  the  linear 


combination  of  Gaussian  orbitals  method  which -bae^ just  -been  developed. 

A  completely  general  form  of  the  spatial  dependence  of  the  potential 

t.  •'  - 

-has  been  used  to  describe  accurately  the  bonding  character  in  the 
tetrahedral  environment.  The  ground  state  charge  densities  we  have1 
determined  agree  very  well  with  experiment  (where  available)  as  do 
the  locations  of  the  valence  bands,  energies,  and  effective  masses. 

A  striking  result  is  that  the  optical  band  gaps  are  underestimated  by 
approximately  30%  or  more,  although  the  general  conduction-band 
topology  is  good.  This  is  an  inherent  difficulty  of  using  the  local 
density  approximation  for  the  exchange  correlation  potential  (a  theory 
of  the  ground  state  properties)  to  describe  the  excitation  energies. 

>—  The  interband  optical  properties  of  Si,  Ge,  GaP,  GaAs,  ZnS,  and 

u'"-" 

ZnSe  -have  been  calculated  using  our  self-consistent  energies  and  wave 

ijj  functions.  Qualitatively  good  agreement  with  experiment  is  found. 

-J 

Agreement  with  experiment  with  regard  to  line  shapes  and  peak  positions 
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can  be  improved  using  an  empirical  energy  dependent  self-energy 
correction  as  appears  in  the  Sham-Kohn  local  density  theory  of 
excitation.  Our  results  clearly  pointed  up  to  the  urgent  need  of  an 
improved  theory  of  the  excitation  energies  in  semiconductors  and 
insulators  which  is  being  actively  pursued, under  the  current  contract 
So  far,  our  work  has  resulted  in  four  extensive  manuscripts,  one 
invited  paper,  and  three  contributed  papers  which  are  enclosed. 
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Progress/Technical  Report 
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by 

C.  P.  Wang 

The  theory  and  computer  programs  have  been  developed  to  calculate 
ab  initio  energies  and  wavefunctions  of  solids  where  the  bonding  has 
important  directional  characteristics.  The  self-consistent  linear 
combination  of  Gaussian  orbitals  method  is  used  with  the  exchange  and 
correlation  potential  included  within  the  local  density  functional 
formulism.  The  crystal  potential  is  expressed  as  a  superposition 
of  spherical  Gaussian  orbitals  centered  at  each  atomic  site,  plus  a 
Fourier  series  expansion  to  accurately  describe  the  bonding  character 
in  the  tetrahedral  environment.  The  effect  of  self-consistency  on 
the  energy  band  gaps  of  diamond  and  zinc  blende  semiconductors  (e.g. 

Si  and  GaAs)  is  currently  under  investigation.  The  application  of 
this  method  to  more  complex  crystal  structures  such  as  Laves  phase 
ZrZn^  (6  atoms/unit  cell)  is  also  in  progress.  Results  will  be  reported 
shortly  at  the  1980  APS  March  meeting. 

For  the  remaining  time  of  the  contract  major  emphasis  will  be  placed 
on  (1)  Application  of  the  method  to  a  series  of  diamond  and  zinc  blende 
semiconductors  in  order  to  make  an  overall  evaluation  of  the  linear  com¬ 
bination  of  Gaussian  orbital  method,  and  (2)  Develop  theory  and  computa¬ 
tional  programs  necessary  to  calculate  the  optical  properties  and  charge 
densities  of  these  materials  in  order  to  make  direct  comparison  of  the 
theoretical  energies  and  wavefunctions  with  the  experimental  measurements. 
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INTRODUCTION 

In  this  paper  we  present  theoretical  results  of  electronic  charge  densities  and  x-ray 
structure  factors  of  diamond  structure  Si  and  Ge,  and  zinc-blende  structure  GaP,  GaAs, 
ZnS  and  ZnSe.  The  calculated  charge  densities  result  from  self-consistent  energy  band 
calculations  for  these  crystals  using  the  method  of  linear  combination  of  Gaussian  orbi¬ 
tals  (LCGO),  with  a  local  density  form  of  exchange-correlation  potential.  The  calcula¬ 
tions  are  fully  ab  initio  with  the  only  input  being  assumed  crystal  structures  and  lattice 
constants.  Further  details  regarding  the  method  may  be  found  in  the  next  section;  and 
the  interested  reader  can  find  a  full  discussion  and  results  for  the  energy  bands,  densi¬ 
ties  of  states,  effective  masses  and  optical  properties  in  Refs.  1  and  2. 

Our  results  discussed  in  Section  III  are  in  good  agreement  with  the  experimental 
structure  factor  measurements  and  indicate  that  local  density  theory  gives  a  good  quan¬ 
titative  representation  of  the  crystal  charge  density  for  these,  and  presumably  other 
materials.  This  is  to  be  expected  as  the  ground  state  charge  density  in  local  density 
theory  results  from  a  variational  principle,  with  the  major  uncertainty  being  the  choice 
of  the  exchange-correlation  functional.  It  appears  that  first-principles  calculations  are 
capable  of  yielding  results  of  sufficient  accuracy  so  that  reliable  interpretations  of  chem¬ 
ical  bonding  in  materials  may  be  made.  Where  theory  and  experiment  deviate,  princi- 
ply  in  the  bond  density  strength  and  asphericity.  improvements  may  be  possible  by 
using  more  realistic  (non-local)  exchange-correlation  functionals. 


METHODOLOGY 


An  effective  one-particle  equation  for  wave  vector  k  and  band  index  n  can  be 
obtained  within  the  Hohenberg-Kohn-Sham  local  density  functional  formalism: 
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Here  the  first  term  is  the  kinetic  energy,  the  second  the  electron-nuclei  Coulomb 
interaction  (nuclei  charge  ZT),  the  third  the  electron-electron  Coulomb  potential,  and 
the  fourth,  the  local  density  exchange-correlation- potential.  The  charge  density  p  (7)  is 
related  to  the  occupied  one-particle  wave  function  0„(k,  7)  by 
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P  (7)  -  ]£  7)0  „  (fc,  7)  (2) 
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which  in  turn  determines  the  electron-electron  Coulomb  and  exchange-correlation 
potentials.  Thus  Eqs.  (1)  and  (2)  need  to  be  solved  self-consistently.  We  have  used 
the  self-consistent  linear  combination  of  Gaussian  orbitals  (LCGO)  method3" 5  which 
has  been  applied  successfully  to  study  the  electronic  properties  of  a  wide  range  of  sim¬ 
ple4-7  and  transition  metals,*-10  and  a  few  covalently  bonded  materials.11,12  Here  we 
will  give  a  brief  discussion  of  the  procedure  that  we  have  adopted  to  calculate  the  self- 
consistent  charge  density  of  representative  semiconductors. 


Linear  Combination  of  Gaussian  Orbitals  Basis 

As  in  the  usual  variational  treatment,  the  crystal  wave  functions  0  „  (k,  7)  are 
expanded  in  terms  of  a  Bloch  basis  set  <f>  iT  (k,  7) 

0„(k.  7)  -  £  C/T.„(k)0/T(k,  7)  (3) 


where  the  standard  LCAO  form  was  used  for  the  basis  function  i  centered  at  the  sublat¬ 
tice  site  r , 
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An  optimum  choice  of  a  basis  set  [«,rC7)]  which  provides  maximum  variational 
freedom  yet  still  maintains  a  manageable  size  of  the  secular  equation  has  been  the  sub¬ 
ject  of  numerous  discussions.13  In  our  study  for  covalently  bonded  diamond  and  zinc- 
blende  materials  we  have  chosen  a  linear  combination  of  atomic  orbitals  minimum  basis 
set  with  an  additional  shell  of  s,  p,  and  d  virtual  Gaussian  type  orbitals  (GTO)  for 
added  variational  flexibility.  The  atomic  orbitals  were  solutions  to  the  atomic  secular 
equation  constructed  with  the  self-consistent  Herman-Skillman  atomic  potential14  and  a 
basis  set  of  GTO.  The  Gaussian  exponents  were  varied  non-linearly  to  minimize  the 
atomic  energies  as  follows:  For  each  orbital  quantum  number  we  choose  a  set  of  even- 
tempered  Gaussians  that  satisfy 
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where  the  two  most  localized  Gaussian  exponents  ai  and  a2  were  varied  non-linearly  to 
minimize  the  corresponding  lowest  state  atomic  energy.  Typically  8,  5,  and  4  GTO 
were  used  to  describe  the  atomic  Is,  2p  and  3d  states,  respectively.  Two  more  extend¬ 
ed  GTO  were  added  with  their  exponents  varied  to  minimize  the  energies  of  the  next 
lowest  states.  This  procedure  was  continued  until  ait  valence  states  of  the  atom  are 
completed.  The  atomic  orbitals  were  then  augmented  by  an  additional  s,  a  and  d  shell 
of  independent  GTO  for  the  variational  freedom  needed  to  describe  the  wave  functions 
in  a  solid.  The  presence  of  these  diffuse  GTO  often  leads  to  approximate  linear  depen¬ 
dence  of  the  basis  set,  hence  the  overlap  matrices  must  be  checked  carefully  against 
negative  or  unphysically  small  eigenvalues  throughout  the  Brillouin  zone.  We  found  it 
helpful  to  truncate  the  tail  of  the  highest  valence  atomic  orbitals  by  setting  minimum 
values  for  the  corresponding  Gaussian  exponents.  Variationally  this  has  little  effect  be¬ 
cause  the  long  tails  of  the  atomic  orbitalrare  strongly  modified  in  a  solid  by  the  overlap 
of  the  wave  functions  on  the  neighboring  atomic  sites.  These  changes  can  only  be 
fitted  with  the  aid  of  the  additional  virtual  orbitals.  The  resulting  basis  set  is  not  com¬ 
plete,  but  does  overlap  both  the  bound  and  low  excited  subspaces  of  the  wave  func¬ 
tions.  It  leads  to  a  dimension  of  36  x  36  for  the  secular  equation  for  Si,  45  x  45  for 
GaP  and  ZnS,  and  54x  54  for  Ge,  GaAs,  and  ZnSe. 


Crystal  Potential 

The  initial  Coulomb  potential  is  constructed  as  a  superposition  of  overlapping 
spherically  symmetric  atomic  potentials.  The  corresponding  overlapping  atomic  charge 
density  is  used  to  evaluate  the  local  density  exchange-correlation  potential  which  we 
chose  to  be  a  Wigner  interpolation  formula15  of  the  form  (energies  are  in  Ry): 

^(7)  -  —  2(3/ir)l/3pl/,Jl.O  +  '  (6) 

where  the  charge  density  p  is  evaluated  at  the  point  7 


Due  to  the  analytic  properties  of  the  Gaussian  orbitals  it  is  customary  to  expand 
the  crystal  potential  in  a  set  of  either  (1)  symmetrized  plane  waves  (SPW) 

K(7)-£KS  £  eM 'r7~7»\  (7) 
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where^the  sum  runs  over  all  group  operations  [R  Id#}  and  stars  of  reciprocal  lattice  vec¬ 
tors  Ks\  or  (2)  overlapping  GTO, 

7(7)-  £  y,Te~m‘' (8) 


centered  at  each  atomic  site. 

The  plane  wave  basis  has  the  advantage  of  being  completely  general,  but  the 
disadvantage  of  slowly  converging  near  the  nuclei  due  to  the  Coulomb  singularity.  For 
example,  even  with  the  aid  of  an  Ewald  type  procedure,  it  is  necessary  to  include  4000 
stars  of  reciprocal  lattice  vectors  to  describe  the  crystal  potential  for  Nis  and  Fe.’  For  a 
more  complex  structure  with  less  symmetry,  such  a  procedure  can  be  prohibitively 
expensive.  The  GTO  on  the  other  hand,  converge  rather  quickly  (less  than  20  GTO 
are  needed  to  describe  an  atomic  potential  to  good  accuracy).  However,  overlapping 
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spherical  GTO  are  too  spherical  around  each  atom  to  describe  the  directional  bond  in  a 
covalent  material.  Although  results  can  be  improved  by  including  higher  angular 
momentum  terms  around  each  atom  and/or  additional  GTO  centered  at  the  interstitial 
tetrahedral  sites,  such  a  procedure  is  somewhat  arbitrary  compared  to  the  Fourier  series 
expansion  where  the  reciprocal  lattice  vectors  are  rigorously  defined.  Thus  we  chose  to 
expand  our  potential  in  a  mixed  basis  of  GTO  and  SPW  as  proposed  by  Euwema.16  The 
self-consistent  Herman-Skillman  neutral  atomic  Coulomb  and  exchange-correlation 
potentials  were  each  fitted  with  18  even-tempered  GTO.  The  four  longest  range  GTO 
which  can  be  represented  by  rapidly  converging  Fourier  series  were  deleted  to  avoid 
excessive  lattice  sums  in  constructing  the  overlapping  crystal  potential.  The  difference 
between  the  exact  atomic  Coulomb  potential  and  the  truncated  GTO  series  was  tabu¬ 
lated  over  a  logarithmic  radial  mesh  and  their  Fourier  coefficients  calculated  numeri¬ 
cally.  In  general  the  crystal  F^(7)  cannot  be.  expressed  as  a  superposition  of  atomic 
Fw(7).  However,  subtracting  out  the  contribution  from  the  overlaping  atomic  GTO 
series  helps  to  eliminate  the  cusp  behavior  near  the  nuclei  to  yield  a  rapidly  converging 
Fourier  expansion  for  the  remainder.  Typically  25  independent  Fourier  coefficients 
were  evaluated  via  a  three  dimensional  least-square  fitting  procedure  based  on  a  sam¬ 
pling  of  400  random  points  in  the  unit  cell.  Having  defined  the  initial  potential  and  the 
basis  functions,  it  is  straightforward  to  evaluate  the  usual  linear  variational  secular 
equation: 

I  \Hir  )a  (k)  -  SITJV  (k)En(k))Cja.n(k)  -  0.  (9) 
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When  expanded  in  terms  of  the  basis  orbitals  ultCr),  the  Hamiltonian  and  overlap 
matrix  elements  are 

Hh  ]a  (?)  -  I  <  ulr (7-  - 7) \H\uj„(7-7tr)>  (10) 
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The  above  two-and  three-center  integrals  for  a  Gaussian  type  basis  can  be  expressed  in 
closed  form.  Results  have  been  published  elsewhere.5  17  The  secular  equation  was 
solved  at  a  set  of  k  points  in  Brillouin  zone  and  a  new  charge  density  was  calculated 
from  Eq.  (2). 

Self-consistency 

The  major  modification  to  the  charge  density  due  to  self-consistency  (charge 
transfer,  change  in  valence  electron  s  and  p  occupation,  etc.)  are  expected  to  occur  in 
the  interstitial  region,  particularly  along  the  tetrahedral  bond.  Therefore  the  Gaussian 
expansion  describing  the  atomic  potential  near  the  nucleus  was  kept  frozen  at  their 
starting  values;  only  the  Fourier  coefficients  were  varied  following  the  procedure  of  Cal¬ 
laway  and  Fry.1®  At  each  iteration,  the  Fourier  coefficients  of  the  input  Coulomb 
potential  VC(K,)  were  evaluated  from  the  output  charge  density  p(K)  of  the  previous 
iteration,  with 

*;(*,>- (12) 
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and  p(Ks)  calculated  analytically  from  the  wave  functions: 

P<*,>-  I  C;TJk)Cjv  l,(k)SlTja(k.-Ks). 
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Here  the  sum  runs  over  all  occupied  bands  (index  n)  and  wave  vectors  k,  CiT  n(k) are 
the  eigenvectors,  and  the  generalized  overlap  matrices 

SltJ<r(k,  K,)  -  / T)e‘K‘  ~<t>jaOc,  T)dir  (14) 

are  simply_the  matrix  elements  of  the  plane  waves  calculated  between  the  basis  func¬ 
tions  <t>lr(k,7).  Ten  special  /c-points19  in  the  1/48'th  of  the  irreducible  Brillouin  zone 
were  used  in  our  iterations  for  a  self-consistent  potential.  In  order_to  speed  up  the  con¬ 
vergence,  p(K,)  was  relaxed  by  mixing  with  50%  of  the  input  p(Ks)  from  the  previous 
iteration  before  it  was  substituted  into  Eq.  (2).  -Once  p( Ks)  has  been  determined  it  is 
straightforward  to  calculate  p(7)  and  hence  the  exchange-correlation  potential  at  the 
selected  points  in  the  unit  cell.  The  contribution  from_  the  overlapping  GTO  are 
subtracted-out  before  the  new  Fourier  coefficients  V^K,)  are  determined  by  a  3- 
dimensional  least-square  fitting  procedure.  Finally  the  plane  wave  contributions  of  both 
the  Coulomb  and  the  exchange-correlation  potential  to  the  new  Hamiltonian  were  cal¬ 
culated  via 

HiT.j„(.k)  -  £(*;(£,)  +  Mtf,))Vj*(*.  *,).  (15) 

where  Slr  Ja  (k,  K, ),  which  needs  to  be  calculated_on!y  once,  is  the  same  matrix  used  in 
Eq.  (13)  to  evaluate  p  ( AT,) .  We  found  that  p  (K,)  converges  to  within  1(T4  a.u.,  and 
the  energies  to  0.02  eV,  at  the  end  of  5  iterations. 

In  the  following  section  we  compare  our  self-consistent  p(K, )  (Eq.  (13))  with 
experimentally  measured  structure  factors.  The  valence  electron  charge  density  maps 
also  shown  were  calculated  with  a  mixed  basis  of  GTO  and  SPW.  The  GTO  were 
obtained  by  fitting  the  atomic  valence  electron  charge  density  p  (7),  while  the  effect  of 
self-consistency  was  included  through  Eq.  (13).  The  sum  over  bands  was  limited  to  the 
valence  electrons  only.  Deformation  densities,  8 p(7),  were  calculated  by  subtracting  a 
superposition  of  spherical  neutral-atom  valence  densities  used  as  starting  values  in  our 
SC  calculations  from  the  final  SC  valence  charge  densities. 


RESULTS  AND  DISCUSSION 

In  this  section  we  present  our  results  for  the  valence  and  deformation  charge  den¬ 
sities,  and  the  x-ray  structure  factors  calculated  according  to  the  procedures  described  in 
the  preceding  section.  All  of  the  experimental  structure  factors  used  for  comparison 
were  reduced  to  zero  temperature  and  corrected  for  anomalous  dispersion  in  order  to 
facilitate  a  meaningful  correlation  between  theory  and  experiment. 

Si  and  Ge 

In  Figs.  1-3  we  show  our  calculated  valence  and  deformation  density  maps  and 
bond  axis  densities,  respectively,  for  Si  and  Ge.  Tables  I  and  II  compares  our  calcu¬ 
lated  structure  factors,  F  Utkl),  with  experiment  and  other  ab  initio  theoretical  results. 
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Fig.  1.  Calculated  self-consistent  valence  charge  density  in  a  portion  of  the  (110)  crystal 
plane  for:  a)  Si,  and  b)  Ge.  The  contours  are  in  units  of  electrons/unit  cell,  and  the 
contour  interval  is  2  electrons/unit  cell. 
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Fig.  2.  Calculated  valence  deformation  charge  density  in  a  portion  of  the  (110)  crystal 
plane  for:  a)  Si.  and  b)  Ge.  The  contours  are  in  units  of  electrons/unit  ceil,  and  the 
contour  intervals  is  1  electron/unit  cell. 
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Fig.  3.  Calculated  p(r)  (solid  lines)  and  8 p(r)  (chain-dotted  lines)  along  the  (111) 
bonding  direction  for:  a)  Si,  and  b)  Ge.  p{r)  is  the  self-consistent  valence  charge  den¬ 
sity,  8 p(r)  is  the  valence  deformation  charge  density  (self-consistent  minus  starting 
overlapping  atomic  charge  densities),  and  r  is  the  fraction  of  the  nearest  neighbor 
(bond)  distance. 


Table  I.  Comparison  of  Theoretical  and  Experimental 
X-ray  Structure  Factors  F(hkl)  for  Si 
(in  Units  of  Electrons/Unit  Cell) 


hkl 

Expt." 

Present 

Results 

Theory* 

Theoryc 

Theory" 

000 

(28.0) 

(28.0) 

(28.0) 

(28.0) 

(28.0) 

111 

15.19 

15.11 

15.07 

15.12 

15.13 

220 

17.30 

17.26 

17.31 

17.28 

17.14 

311 

11.35 

11.37 

11.41 

11.33 

11.02 

222 

0.38 

0.25 

— 

0.34 

0.38 

400 

14.89 

14.92 

14.92 

14.88 

14.70 

331 

10.25 

10.17 

10.15 

10.20 

9.94 

422 

13.42 

13.37 

13.37 

13.36 

13.30 

333 

9.09 

9.07 

9.07 

9.02 

8.92 

511 

9.11 

9.08 

9.08 

9.08 

8.98 

440 

12.08 

12.04 

- 

12.04 

12.00 

“Experiment:  Aldred  and  Hart.  Ref.  27.,  except  for  1222)  which  is  from 
Ref.  28. 

*LCAO  theory:  Heaton  and  Lafon.  Ref.  12. 

‘SCOPW  theory:  Raccah.  et  ai.  Ref.  29 
"Pseudopotentiai  theory:  Zunger  and  Cohen.  Ref.  21. 


ECTRONS/UNIT  CELL) 


Table  II.  Comparison  of  Theoretical  and 
Experimental  X-ray  Structure 
Factors  F(hkl)  for  Ge 
(in  Units  of  Electrons/Unit  Cell) 


hkl 

Expt.a 

Present 

Results 

SCOPW4 

000 

(62.0) 

(62.0) 

(62.0) 

111 

39.42 

38.83 

38.95 

220 

47.44 

47.23 

47.26 

311 

31.37 

31.29 

31.21 

222 

0.26 

0.22 

0.48 

400 

40.50 

40,56.. 

40.54 

331 

27.72 

27.39 

27.53 

422 

36.10 

35.91 

36.00 

333 

24.50 

24.35 

24.34 

511 

— 

24.35 

24.38 

440 

32.34 

32.23 

32.23 

“Matsushita  and  Kohra,  Ref.  24. 

^Raccah.  et  at.,  Ref.  29. 

Figures  1-3  clearly  show  the  covalent  nature  of  the  bonding  in  Si  and  Ge,  with  the 
characteristic  charge  pile-up  near  the  bond  center  for  both  materials.  It  is  further 
interesting  to  note  from  the  deformation  densities  shown  in  Figs.  2  and  3  how  charge  is 
pulled  into  the  bond  compared  with  the  situation  for  overlapping  charge  densities— 
much  of  the  covalency  results  from  charge  redistribution,  not  merely  from  wave  func¬ 
tion  overlap. 

The  valence  and  deformation  maps  are  in  good  agreement  with  previous  ab  initio 
calculations20, 12  and  the  experimental  densities  constructed  by  Yang  and  Coppens22  for 
Si.  The  characteristic  features  of  the  valence  density  being  elongated  along  the  bond, 
and  opposite  for  the  deformation  density,  are  reproduced  in  the  calculated  maps.  How¬ 
ever.  Yang  and  Coppens22  find  a  valence  bond  maximum  of  27.6,  while  we  find  21.3 
(electrons/unit  cell)  for  Si.  Our  value  is  very  close  to  the  theoretical  value  found  by 
Hamann20  (  22.2  electrons/unit  ceil)  but  somewhat  smaller  than  that  found  by  Zunger 
and  Cohen21  (24.0  electrons/ unit  cell).  Recently,  Scheringer22  has  pointed  out  that  the 
Yang  and  Coppens22  value  may  be  overestimated  by  as  much  as  3.6  electrons/ unit  cell 
due  to  inaccuracies  in  four  of  the  measured  high  momentum  structure  factors.  For  Ge, 
where  experimental  structure  factor  data  exists  for  only  a  relatively  few  number  of 
reflections24  2-  (sin  0/A  <  0.50  /l-1).  quantitative  results  for  the  bonding  have  been 
limited.26 

Calculated  and  measured  structure  factors  for  Si  and  Ge  are  compared  in  Tables  I 
and  II,  respectively.  For  silicon,  except  for  the  (222)  forbidden  reflection,  the  experi¬ 
mental  results  of  Aldred  and  Hart27  who  used  the  Pendellosung  fringe  method  are 
shown  in  Table  I.  For  the  (222)  reflection  in  Si  the  value  shown  was  chosen  from  the 
work  of  Roberto,  Batterman  and  Keating,28  corrected  to  zero  temperature  (see  Ref. 
28).  The  experimental  results  of  Matsushita  and  Kohra24  for  Ge  were  chosen  for  com¬ 
parison  in  Table  II.  These  measurements  were  made  using  the  half-width  of  the 
Bragg-case  diffraction  curves  measured  in  the  triple  crystal  arrangement  using  Cu-Ka 


radiation.  Tables  1  and  II  also  show  the  present  theoretical  structure  factor  results, 
together  with  calculated  values  of:  Heaton  and  Lafon12  (for  Si)  who  used  a  self- 
consistent  LCAO  method  similar  to  ours  but  with  a  Kohn-Sham  (a  -2/3)  exchange- 
correlation  potential  and  a  smaller  basis  set;  Raccah,  et  al.29  (for  Si  and  Ge)  who  used  a 
self-consistent  orthogonalized  plane  wave  (SCOPW)  approach  also  with  a  »■  2/3;  and 
Zunger  and  Cohen21  (for  Si)  who  used  the  first-principles  self-consistent  pseudopoten¬ 
tial  method  with  the  Hedin-Lundqvist30  form  of  local  density  exchange-correlation 
potential. 

Our  structure  factor  results  for  Si  compare  extremely  well  with  the  measurements 
of  Aldred  and  Hart.27  Most  of  the  theoretical  results  also  agree  very  well,  especially  the 
present  results  and  the  other  SC-LCAO12  and  the  SCOPW29  calculations.  There  are 
more  substantial  disagreements  between  the  SC-pseudopotential21  results  and  the  other 
calculations,  and  with  the  measurements.  The  exception  is  for  the  (222)  reflection 
where  the  present  results  are  considerably  smaller  than  experiment,  while  the  SC-OPW 
and  SC-pseudopotential  results  are  much  closer.  Since  the  calculations  of  Heaton  and 
Lafon12  and  our  own  were  done  in  a  similar  manner  except  for  the  choice  of  exchange- 
correlation  potential  and  basis  functions,  the  very  close  agreement  between  the  two  sets 
of  results  seems  to  indicate  that  (at  least  for  Si)  the  self-consistent  charge  density  is 
rather  insensitive  to  these  choices. 

Our  structure  factor  results  for  Ge  also  compare  very  well  with  the  measurements 
of  Matsushita  and  Kohra24  and  with  the  SCOPW29  calculations.  Here  too  our  calculated 
value  of  F(222)  is  smaller  than  the  measured  value,  though  closer  than  was  the  case  for 
Si.  The  SCOPW  value  of  F(222)  for  Ge  is  almost  a  factor  of  two  larger  than  experi¬ 
ment. 

GaP.  GaAs,  ZnS  and  ZnSe 

Figures  4-9  show  our  calculated  valence  and  deformation  maps  and  bond  axis 
densities  for  the  zinc-blende  compounds.  Tables  III-V  compare  the  absolute  values  of 
our  calculated  structure  factors  with  experiment  and  with  other  ab  initio  results,  while  in 
Table  VI  we  present  our  results  for  both  the  real  and  imaginary  parts  of  F(  hkl)  for  all 
four  compounds.  It  should  be  noted  that  the  density  maps  for  ZnS  and  ZnSe  include 
the  3d-states  of  Zn  and  Se,  since  the  Zn  3d-states  fall  in  the  upper  valence  band  region 
(see  Ref.  1),  while  for  the  gallium  compounds,  the  much  lower-lying  Ga  and  As  3d- 
states  were  kept  in  the  core.  Density  maps  with  the  Zn  and  Se  3d-states  kept  in  the 
core  may  be  found  in  Ref.  1. 

The  valence  density  maps  for  ail  four  materials  in  Figs.  4  and  5  show  a  charge 
build-up  towards  the  anion  sites  (P,  As,  S,  or  Se).  However,  the  deformation  maps 
and  bond  axis  densities  shown  in  Figs.  6-9  are  even  more  instructive.  The  deformation 
maps  clearly  show  that  the  bonding  induced  by  the  crystalline  environment  is  both  ionic 
and  covalent.  In  fact,  the  deformation  maps  are  qualitatively  similar  to  Si  and  Ge 
insofar  as  they  show  a  bond  charge  buildup  elongated  perpendicular  to  the  bond  axis, 
but  with  the  important  difference  that  the  deformation  bonding  charge  is  shifted 
towards  the  anion  site.  From  Figs.  8  and  9  it  can  be  seen  that  the  bond  axis  deforma¬ 
tion  charge  densities  peak  towards  the  anion  site,  but  the  bond  axis  valence  densities 
are  shifted  even  further  towards  the  anion,  indicating  that  a  large  part  of  the  valence 
charge  asymmetry  is  due  to  atomic  wave  function  overlap  effects. 


Table  III.  Comparison  of  Theoretical 
and  Experimental  Absolute  Values 
of  The  X-ray  Structure  Factors 
\F(hkl)\  for  GaAs 
(in  Units  of  Electrons/Unit  Cell) 


hkl  Expt. 


Present 

Results 


"Experimental  results  of  Matsushita  and 
Havashi  (Ref.  31)  corrected  for  anomalous 
dispersion  (/'-  -  1.31)  and  temperature 
<fl  -  0.595 

^Experimental  results  of  Matsushita  and 
Havashi  <  Ref  31)  corrected  for  anomalous 
dispersion  -  1.01)  and  temperature 

IS  -0.629 


Table  IV.  Comparison  of 
Theoretical  and  Experimental 
Absolute  Values  of  The 
X-ray  Structure  Factors 
\F(hkl)  \  for  GaP 
(in  Units  of  Electrons/Unit  Cell) 


hkl  Expt."  Prese,m 
Results 


2 


Table  V.  Comparison  of  Theoretical  and 
Experimental  Absolute  Values  of  the 
X-ray  Structure  Factors  \F(hkl)\  for  ZnSe 
(in  Units  of  Electrons/ Unit  Cell) 


“Ref.  29. 


Table  VI.  Calculated  Real  and  Imaginary  Parts  of  the  X-ray  -  • 

Structure  Factors  F(hkl)  for  the  Zinc-Blende 
Compounds.  The  origin  is  at  a  cation  site,  and  the 

units  are  electrons/unit  cell.  l 


In  Table  III  we  compare  our  calculated  values  of  \F{hkl)\  for  GaAs  with  the 
experimental  results  of  Matsushita  and  Hayashi31  corrected  for  thermal  vibrations  and 
anomalous  dispersion.  The  experimental  data  were  determined  using  the  half-width  of 
the  Bragg-case  diffraction  curves  measured  in  the  triple  crystal  arrangement  using  Cu- 
Ka  radiation  (similar  to  Matsushita  and  Hayashi's  experiments28  on  Ge  described 
above).  For  the  case  of  polyatomic  materials  it  is  considerably  more  uncertain  how  to 
correct  for  anomalous  dispersion  (fj)  and  determine  the  Debye-Waller  correction 
(e~u).  Different  values  of  fj  exist  in  the  literature  for  Ga  and  As,  and  the  value  of  B 
(M  —  5(sin0/X)2)  depends  on  the  choice  of  fj.  Matsushita  and  Hayashi3'  discuss 
their  results  in  terms  of  two  sets  of  fj,  from  Cromer32  and  from  Cromer  and  Liber¬ 
man,33  respectively,  without  making  a  clear-cut  choice  J>etween  the  sets.  We  have 
therefore  reduced  their  experimental  data  by  considering  f  as  the  Ga  and  As  mean  for 
both  sets:  a)  /*  —  -1.31,  B  *  0.595  A1,  and  b)  f  —  -1.01,  B  *■  0.629  A2.  Both  sets  of 
reduced  experimental  data  are  shown_in  Table  III  and  compared  with  our  theoretical 
results.  Given  the  uncertainty  in  the  f  and  B  corrections  to  the  experimental  data,  and 
the  quoted  —1%  accuracy  in  the  measurements,  the  agreement  between  theory  and 
experiment  is  remarkably  good. 

lino  and  Ishigaki34  have  measured  the  structure  factors  of  GaP  using  two  different 
methods:  a)  the  angle  dispersive  method  using  monochromatized  Cu-Ka  x-rays  and,  b) 
the  white  x-ray  diffraction  method,  both  performed  on  powder  samples.  The  angle 
dispersive  method  results,  being  more  accurate  (see  Ref.  34),  are  compared  with  our 
calculations  in  Table  IV.  The  experimental  results  shown  have  been  corrected  to  zero 
temperature  using  a  value  of  B  —  0.84  obtained  from  Fig.  3  of  Ref.  34 
[a  log  ( | F ;Xpl/ 1 F iheorv I )  plot  using  Cromer’s32  values  for  f'Ga  and  fj> ].  As  for  the  GaAs 
case,  there  is  some  uncertainty  in  the  f  and  B  values,  and  the  quoted  experimental 
accuracy  is  --1-3%.  Within  these  experimental  uncertainties  our  results  for  GaP  again 
agree  very  well  with  the  measurements. 

For  ZnSe  we  compare  our  results  with  the  experiments  and  SCOPW  calculations 
of  \F(hkl)\  of  Raccah,  et  al.  in  Table  V.  The  measurements  were  done  on  a  powder 
sample  with  the  thermal  and  anomalous  diffraction  corrections  taken  into  account  by 
Raccah,  et  al.19  (they  present  the  corrected  results),  with  a  quoted  experimental  uncer¬ 
tainty  of  several  percent.  Our  calculations  agree  quite  well  with  the  measurements, 
with  the  most  substantial  disagreement  being  for  |F(222)|.  Our  results  also  agree  very 
well  with  those  resulting  from  the  SCOPW  calculations. 


CONCLUSIONS 

We  can  conclude  from  the  theoretical-experimental  comparisons  made  in  this 
paper  that  local  density  theory  gives  a  good  description  of  charge  densities  in  semicon¬ 
ductor  systems.  The  major  discrepancy  between  theory  and  experiment  appears  to  be 
that  the  theory  underestimates  the  charge  density  asphericity,  as  evidenced  by  the  too 
small  vlaues  of  F(222)  and  bond  charge  for  Si  and  Ge.  This  seems  to  be  a  general 
result  of  many  SC  calculations  for  charge  and  spin  densities  and  is  due  to  the  underly¬ 
ing  local  density  approximation.  Some  examples  for  metals  are  shown  in  Table  VII. 
Improvements  in  the  exchange-correlation  functional  (e.g.,  non-local  corrections)  may 
resolve  these  subtle  but  interesting  disagreements  between  theory  and  experiment. 


Table  VII.  Directional  Anisotropy  of  Charge  and  Spin  Densities 
for  Transition  Metals  by  Comparing  the  Ratios  of  Structure 
Factors  for  Wave  Vectors  of  the  Same  Magnitude. 


Element 

Structure 

Factors 

Wave 

Vectors 

Theory 

Experiment 

Ni 

Spin 

(3.3,3) 

(5.1.1) 

2.179° 

3.028* 

Fe 

Spin 

(3,3.0) 

(4.1.1) 

2.35 1 c 

2.848* 

Fe 

Charge 

(3,3,0) 

(4.1.1) 

1.0025' 

1.010' 

1.01 1' 

V 

Charge 

(3.3,0) 

(4.1.1) 

1.0039* 

1.024* 

1.0085' 

"Ref  8.  *Ref  35.  'Ref  9,  “'Ref.  36. 
'Ref.  37,  rRef  38,  *Ref  10.  *Ref.  39. 
'Ref.  40. 


However,  it  is  clear  from  the  present  study  that  good  quality  ab  initio  calculations  are 
capable  of  giving  useful  and  reliable  results  for  the  bonding  and  chemistry  of  covalently 
bonded  semiconductors. 
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